How many diffusion gradient directions are required for HARDI?

Introduction: While the minimum number of directions required for robust diffusion tensor imaging (DTI) has been extensively studied [1-4], the equivalent studies have not been performed for high angular resolution diffusion-weighted imaging (HARDI) methods. This is due in part to the fact that many HARDI methods are available, each with their own requirements and parameters, making it difficult to select an appropriate outcome measure. It is also difficult to replicate the data acquisition strategies that have been employed in the original DTI studies, which consisted of acquiring a number of schemes with different numbers of directions and repeats, while keeping the total number of images constant (e.g. 10×6, 5×12, 3×20, 2×30, 1×60). This strategy was possible with DTI due to the known theoretical minimum number required (i.e. 6), and the fact that this number was sufficiently low to keep the total scan time feasible. In this study, we use an approach based on sampling theory to determine the minimum number of DW directions required for HARDI experiments, independent of the particular algorithm used. Theory: The diffusion-weighted (DW) signal measured as a function of orientation (referred to here as the signal profile) can be viewed as a 2D waveform defined over the sphere. The number of unique sampling directions needs to be sufficiently large to provide a proper characterisation of all the features of the signal profile. This can be appreciated using a Cartesian Fourier analogy: the Nyquist criterion states that the minimum sampling frequency required to fully characterise a waveform is twice the maximum frequency component contained in that waveform. Similarly, the minimum sampling frequency in the angular domain required to fully characterise a waveform over the sphere is determined by its maximum angular frequency content. The spherical harmonic (SH) series provides a rigorous framework to characterise these angular frequencies. If the signal profile can be shown to contain no significant SH terms beyond a given harmonic order lmax, then the minimum number N of unique directions required to fully characterise it is given simply by the number of SH coefficients contained in the SH series up to lmax, N = 1⁄2(lmax+1)(lmax+2). Note that only even harmonic orders are needed due to the symmetry of the diffusion process. The signal profile will be sharpest, and hence contain the highest angular frequency content, in regions containing a single coherently oriented fibre population. This can be more readily appreciated by considering the problem within the spherical deconvolution framework [5]: the signal profile within a crossing fibre region is given by the spherical convolution of the fibre orientation distribution (FOD) with the signal profile for a single fibre orientation. According to the convolution theorem, convolution over one domain (the sphere) is equivalent to multiplication in the other domain (angular frequency). Hence, the signal profile measured in a crossing fibre region cannot contain higher angular frequencies than were already present in the single-fibre signal profile. In addition, the signal profile in single-fibre regions can be assumed to be axially symmetric. This is a substantial advantage, since in this case only one free parameter (the m=0 coefficient) needs to be estimated per harmonic order l instead of 2l+1 for the non-axially symmetric case. Methods: The frequency content of the signal profile was estimated in two ways. First, the noiseless signal profile was simulated for a single fibre population using the diffusion tensor model, assuming a b-value of 3000 s/mm and literature values for the mean ADC and FA of 500 × 10 mm/s and 0.7 respectively [6]. The m=0 (axially symmetric) SH coefficients of the simulated DW signal were then estimated using a least-squares fit. Second, the signal profile was measured directly from a high-quality in-vivo data set, acquired from a healthy volunteer using a 3T Siemens Trio MRI system. Acquisition parameters were: b = 3000 s/mm, TE = 110 ms, 3 mm isotropic voxels, 30 contiguous slices, 500 uniformly distributed DW gradient directions, 25 b=0 images, scan time: ~36 minutes. The analysis was restricted to the set of voxels with FA greater than 0.6, since these were assumed to contain a single fibre population. Within each voxel, the m=0 (axially symmetric) SH coefficients of the signal profile were estimated using a maximum likelihood estimation procedure. The likelihood of the data within each single-fibre voxel, given the set of m=0 SH coefficients {sl} and the noise σ is given by: